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Abstract 

Although unstable, the Si(337) orientation has been known to appear in diverse experimental 
situations such as the nanoscale faceting of Si(112), or in the case of miscutting a Si(113) surface. 
Various models for Si(337) have been proposed over time, which motivates a comprehensive 
study of the structure of this orientation. Such a study is undertaken in this article, where 
we report the results of a genetic algorithm optimization of the Si(337)-(2 x 1) surface. The 
algorithm is coupled with a highly optimized empirical potential for silicon, which is used as an 
efficient way to build a set of possible Si(337) models; these structures are subsequently relaxed 
at the level of ab initio density functional methods. Using this procedure, we retrieve the (337) 
reconstructions proposed in previous works, as well as a number of different ones. 
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1 Introduction 

Silicon surfaces have been at the foundation of the semiconductor industry for several decades. 
With the advent of nanotechnology, the interest in the atomic structure of various semiconductor 
surfaces has expanded because of their potential to serve as natural and inexpensive templates 
for growing various types of nanostructures. Understanding the growth of nanoscale entities on 
various substrates depends, at least in part, on knowing the atomic structure of the substrate. As 
the low-index crystal surfaces are well understood, the high-index orientations are now progressively 
gaining importance as their richer morphological features may lead to more advanced technological 
applications. To date, there are numerous stable high-index orientations reported in the literature; 
in the case of Silicon and Germanium, these stable orientations have been summarized in Ref. [Q. 
The unstable high-index orientations are important in their own right, as they often give rise to 
remarkable periodic grooved morphologies which may be used for the growth of surface nanowires. 
The wonderful morphological and structural diversity of high-index surfaces can be appreciated, for 
instance, from the work of Baski, Erwin and Whitman, who investigated systematically the surface 
orientations between Si(OOl) and Si(lll) 

An interesting case of unstable surface is Si(337), whose atomic structure is the subject of the 
present article. As a historical account, we note that the Si(337) orientation was concluded to 
be stable in early studies by Ranke and Xing Gardeniers et al. 4 , and Hu et al. 5 j. This 
conclusion about the stability of Si(337) was consistent with the reports of Baski and Whitman, who 
showed that Si(112) decomposes into quasi-periodic Si(337) and Si(lll) nanofacets [5], and that 
Si(337) is a low-energy orientation [7]. However, further studies by Baski and Whitman revealed 
that Si(337) itself facets into Si(5 5 12) and Si(lll) demonstrating that Si(337) was, in fact, 
unstable. Since Si(337) is not stable and Si(5 5 12) is [H1E3> the clean Si(112) orientation should 
facet into Si(5 5 12) and Si(lll). However, as explained in Ref. nanofacets are too narrow 
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to form complete Si(5 5 12) units, and mostly settle for the nearby Si(337) orientation. A more 
recent report shows another situation where Si(337) arises rather than the expected Si(5 5 
12). For a Si(113) surface that is miscut towards [111], the large (5 5 12) reconstruction has not 
always been found. Instead, a Si(337) phase has been seen to coexist with Si(113) terraces |llj . 
The high resolution transmission electron microscope (HRTEM) images in Ref. jllj leave no doubt 
that the periodicity of those nanofacets corresponds to the (337) orientation. The reason for the 
quasi-stability of Si(337) reported in earlier works QJ could be the curvature or the size of 
substrates used in those investigations; this explanation is also consistent with Ref. |llj . where the 
issue of sample curvature is acknowledged. 

With one exception [5], atomic-scale models for Si(337) were not sought separately, but rather as 
structural parts of the Si(5 5 12)-(2 x 1) reconstruction H21 H3J EH • As shown in Ref. the Si(5 
5 12) unit cell consists of two Si(337) units and one Si(225) unit. The first model reconstruction 
proposed for Si(5 5 12) [3] appeared somewhat corrugated when viewed along the [110] direction; 
so did a second model by Ranke and Xing |12j . On the other hand, HRTEM measurements |13] 
showed a flatter profile for Si(5 5 12), and different model reconstructions were proposed |13l I14j 
in order to account for the observed flatness. Total-energy density functional calculations for the 
energy of Si(5 5 12) have only recently been published, and suggest that the latest models [T3*l ITi] 
have lower energies than the earlier proposals [HI El • 

The reason structural studies of high-index surfaces can be very lengthy originates from the way 
model reconstructions are proposed, which is based on physical intuition in interpreting scanning 
tunneling microscope (STM) images. As seen above for the case of Si(5 5 12), the intuitive ways 
are not very robust for high- index surfaces, since usually more models are required when further 
experiments are reported. Given the rather large number of low-energy structures that high-index 
surfaces can have (see, e.g., the case of Si(105)[T51 ITB]). heuristic approaches pose the risk of 
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missing the actual physical reconstruction when considering model candidates for comparison with 
the experimental data. To alleviate this risk, we have recently addressed the reconstruction of 
semiconductor surfaces as a problem of stochastic optimization ^H], and developed two different 
global search algorithms for the purpose of determining the most favorable surface structure |161ll7j . 
We have used these global search procedures for studying Si(105) ^3 El and Si(114) ^H]) which 
are both stable orientations for silicon. Although Si(337) is not a stable orientation, its atomic 
structure can still be pursued with global search algorithms because the only input required is the 
surface periodicity, which is available from STM |5] and LEED [3] measurements. 

In this article we report the results of such a global optimization for the structure and energy 
of Si(337) using empirical potential and density functional theory calculations. As we shall see, 
the optimization retrieves most of the previously proposed models jSJ |§J 1131 I14j . as well as a 
number of other reconstructions that could be relevant for experimental situations where Si(337) 
nanofacet arises — other than the two (337) unit cells that appear as part of the (5 5 12) structure. 
The remainder of this paper is organized as follows. In the next section we briefly describe the 
computational methods used in this study. In Section 3, we present our model structures for the 
Si(337) surface classified by the number of atoms in the computational cell and their surface energy. 
In Section 4, we compare our results with other reported Si(337) reconstructions. Our conclusions 
are summarized in the last section. 

2 Methods and Computational details 

The Si(337) simulation cell consists of two (337) bulk truncated unit cells of dimensions a\/8.375 x 
avOJ3xav67, (a is the lattice constant of Si), stacked along the y-direction, [110] . These dimensions 
are calculated from the crystal geometry, and are in agreement with the STM images published in 
Ref. 0, where (337) units are shown as part of a stable reconstructed surface, Si(5 5 12)-(2 x 1). 
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Because the (337) plane has two-dimensional periodicity along two non- orthogonal directions, the 
rectangular unit in Fig. ^ must be subjected to shifted periodic boundary conditions JHJ) with the 
shift in the y-direction equal to a/y8. We build a database of structures sorted out by the number 
of atoms in the simulation cell and by their surface energy calculated with the highly optimized 
empirical potential developed by Lenosky et al. 20 . The number of atoms is varied around the 
initial value of n = 268, in order to exhaust all the distinct possibilities that can occur during the 
search for low energy reconstructions. 

The global optimization involved in building the structure database for Si(337) has been per- 
formed with the genetic algorithm adapted for crystal surfaces |17j . This algorithm is based on 
an evolutionary approach in which the members of a generation (pool of models for the surface) 
evolve with the goal of producing the best specimens, i.e. lowest energy reconstructions. The 
evolution from one generation to the next takes place by mating, which is achieved by subjecting 
two randomly picked structures {parents) from the pool to an operation that combines their sur- 
face features and creates a child configuration. This child structure is relaxed and considered for 
inclusion in the pool, depending on its surface energy. We use two variants of the algorithm, one in 
which the number of atoms is kept constant for all members of the pool at all times (the constant-n 
version), and another one in which this restriction is not enforced (the variable- ra version). The 
implementation of this genetic algorithm has recently been described in various degrees of detail 
in Refs. |17l 1181 l2*T] : for this reason, we will not expand upon it here, but refer the reader to those 
works. Since the empirical potential may not give a very reliable energy ordering for a database of 
structures, we study not only the global minima given by HOEP for different values of n, but also 
most of the local minima that are within 15 meV/A 2 from the lowest energy configurations. After 
creating the set of Si(337) model structures via the genetic algorithm, we recalculate the surface 
energies of these minima at the level of density functional theory (DFT). 
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The DFT calculations where performed with the plane- wave based PWscf package |22j . using 
the Perdew-Zunger exchange-correlation energy. The cutoff for the plane-wave energy was set 
to 12 Ry, and the irreducible Brillouin zone was sampled using 4 /c-points. The equilibrium bulk 
lattice constant was determined to be a =5.4lA, which was used for all the surface calculations in 
this work. The simulation cell has the single-face slab geometry with a total silicon thickness of 
10A , and a vacuum thickness of 12 A. The Si atoms within a 3A-thick region at the bottom of 
the slab are kept fixed in order to simulate the underlying bulk geometry, and the lowest layer is 
passivated with hydrogen. The remaining Si atoms are allowed to relax until the force components 
on any atom become smaller than 0.025 eV/A. 

The surface energy 7 for each reconstruction is determined indirectly, by first considering the 
surface energy 70 of an unrelaxed bulk truncated slab, then by calculating the difference A7 = 7—70 
between the surface energy of the actual reconstruction and the surface energy of a bulk truncated 
slab that has the bottom three layers fixed and hydrogenated. This indirect method for calculating 
the surface energies at the DFT level was outlined, for instance, in Ref. |24| . The energy of the 
bulk truncated surface shown in Fig. ^was found to be 70 = 137.87 meV/A 2 . At the end of the 
genetic algorithm search, we obtain a set of model structures which we sort by the number of atoms 
in the simulation cell and by their surface energy. These model structures are presented next. 

3 Results 

There are four possibilities in terms of the number of atoms in the slab that yield distinct global 
energy minima of the Si(337) periodic cell shown in Fig.^ This has been determined by performing 
constant-n genetic algorithm optimization for computational slabs with consecutive numbers of 
atoms n (264 < n < 272), and identifying a periodic behavior of the lowest surface energy as a 
function of n. This procedure of determining the symmetry distinct numbers of atoms in the slab 
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has been detailed in Ref. |16j . The HOEP global minima, as well as selected local minima of the 
surface energy for different numbers of atoms are summarized in Table HI along with the density of 
dangling bonds per unit area and the surface energy computed at the DFT level. 

As a general comment, the DFT energy ordering does not coincide with that given by HOEP, 
indicating that the transferability of HOEP to Si(337) is not as good as in the case of Si(OOl) 
and Si(105); to cope with this transferability issue, we have considered more local minima (than 
listed in TableQJ when performing DFT relaxations, as mentioned in Section 2. In the terminology 
of potential energy surface (PES) theory, we sample the main basins using the genetic algorithm 
coupled with HOEP, then recalculate the energetic ordering of these basins using DFT. From the 
table it is apparent that the least favorable number of atoms is n = 267 (modulo 4) at both 
the HOEP and DFT levels. Therefore, we focus on describing here the reconstructions that have 
numbers of atoms n = 266, 268, and 269. 

For n = 266 the best model that we obtained at the DFT level is made of a fused assembly of 
dimers, tetramers and slightly puckered honeycombs, followed (in the [776] direction) by a row of 
rebonded atoms. These individual atomic-scale motifs have been reported previously for different 
surfaces (e.g., [51 1141 [2*5] ). and depicted for convenience in Fig. |2J Their complex assembly shown 
in Fig. EKa) has a surface energy of 94.47 meV/A 2 ; the reconstruction Of a) has a corrugated aspect 
when viewed from the [110] direction, which may account for its relatively high surface energy (refer 
to Table P). Higher-energy models with n = 266 are illustrated in Fig. Etb,c,d), and have surface 
energies ranging from 2.75 meV/A 2 to about 7 meV/A 2 above the surface energy of the model|Hta). 
Dimers (D) and tetramers (T) can be clearly identified on model [2Kb), while configuration E[c) is 
characterized by the presence of honeycombs (H) and rebonded atoms (R). In model Etd) the only 
features from the set in Fig. [21 that can be separately identified are the rebonded atoms R. Nearly 
degenerate with structure 3(d), we find two other models [3(e) and Off ), with different tilting of 
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the dimers belonging to tetramer groups] that contain 6- member ring 7r-chains The 6- member 
rings are labelled by 6-r in Fig. Efe,f), and are supported by tetramer- like features (denoted by t 
in Fig. I2J) whose dimers are made of fully coordinated atoms. 

The optimum number of atoms is n = 268, as indicated in Table [Q at both HOEP and DFT 
levels -coincidentally, this is also the number of atoms that corresponds to two complete (337) bulk 
truncated primitive cells. We have further verified this number by performing a variable-n genetic 
algorithm started with all the members in the pool initially having n = 267 atoms. The algorithm 
has found the correct number of atoms (n = 268) and the HOEP global minimum in less than 1000 
operations. While different starting configurations may change this rather fast evolution towards 
the global minimum, our tests indicate that for simulation cells of similar size several thousand 
genetic operations are sufficient when using the implementation described in Ref. ^7]. The lowest 
energy model as given by DFT is shown in Fig. Eta) and consists of dimers D, rebonded atoms 
R and honeycombs H, in this order along [776]. We note that the number of dangling bonds is 
smaller at the DFT level (refer to Table QJ): the reason for the dangling bond reduction is the 
flattening of the honeycombs (not planar at the level of empirical potentials), which creates two 
additional bonds with the subsurface atoms. The genetic algorithm search retrieves another low- 
energy model, in which the dimers are displaced towards and combine with the rebonded atoms, 
forming tetramers T [Fig. Efb)]. The two models in Figs. E[a) and (b) are nearly degenerate, 
with a surface energy difference of ~ 0.5 meV/A 2 relative to one another. A different number 
and ordering of the structural motifs that are present in Fig. Eta) gives rise to a noticeably larger 
surface energy. Specifically, the arrangement of dimers D and honeycombs H [Fig. EI C )] increases 
the surface energy by 4 meV/A 2 relative to the model Ufa). An even higher-energy reconstruction 
made up of tetramers and rebonded atoms is shown in Fig. 0fd). 

When the number of atoms in the simulation cell is n = 269, we find that the ground state is 
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characterized by the presence of pentamers, subsurface interstitials, and rebonded atoms [refer to 
Fig. Eta)]- Interestingly, this structure is nearly flat, with the rebonded atoms being only slightly 
out the plane of the pentamers. The pentamers are "supported" by six-coordinated subsurface 
interstitials. For the Si(337) cell with n = 269 atoms, we find a very strong stabilizing effect of the 
pentamer-interstitial configuration: the next higher-energy structure [Fig. E^b)] is characterized by 
adatoms and rebonded atoms, with a surface energy that is ~5 meV/A 2 higher than that of the 
ground state in Fig. O^a). We note that subsurface interstitials have long been reported in the 
literature to have a stabilizing effect for the reconstruction of Si(113) |24| . Next, we compare the 
models described above with previously published works on the Si(337) structure. 

4 Discussion 

To our knowledge, the first model of the Si(337) surface was proposed by Baski, Erwin and Whitman 
in their original work on the structure of Si(5 5 12) The model contained surface tetramers and 
6-member ring 7r-chains, and was similar to the structures depicted in Fig. Ote,f). We observed a 
strong relaxation of the 7r-bonded chain, in which one of the atoms of each 6-ring protrudes out 
of the surface and settles well-above its neighbors. The surface energy of the model |2tf) is not 
particularly favorable, although it is the lowest of all the models containing 6-ring 7r-chains that 
we have investigated. Heuristically speaking, if the under-coordinated atoms of the 6-rings are 
removed (4 atoms per (2 x 1) unit cell), then one would obtain the model of Hu et al. [S], which is 
made of dimers and tetramers. The genetic algorithm with n = 266 (modulo 4) has also retrieved 
this model, which is shown in Fig. E^b) after DFT relaxation. 

Models|3Ib,e,f) have a high density of dangling bonds (db) per unit area, 10db/a 2 ^/l6?75. While 
the dangling bond density may seem to account for the large surface energy of these models, a close 
inspection of Table [I] shows that, in fact, there is no clear correlation between the dangling bonds 
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and the surface energy. For instance, at the same number of unsatured bonds per area, models 
E[b,e,f) span a 4 meV/A 2 range in surface energy. Furthermore, if we take the DFT model |2Ka) as 
reference, we note that an increase by one db per unit cell [fromEfa) to|Hfc)] results in a surface 
energy increase of ~7 meV/A 2 ; on other hand, an increase of seven db in the simulation cell [from 
Ota) to|3tb)] increases the surface energy by less than 3 meV/A 2 . Similar arguments for the lack 
of correlation between missing bonds and surface energy can be made by analyzing structures with 
other numbers of atoms in Tabled At the optimum number of atoms (n = 268), we find that 
with the same dangling bond density, models 0fa,b,c) span a 4 meV/A 2 surface energy range: this 
range is so large that it includes the global minimum as well as higher surface energy local minima 
that may not be observable in usual experiments. These findings are consistent with the statement 
that the minimization of surface energy is not controlled solely by the reduction of the dangling 
bond density, but also by the amount of surface stress caused in the process Given that the 
heuristic approaches may only account for the dangling bond density when proposing candidate 
reconstructions, the need for robust minimization procedures (e.g., |16l I17j ) becomes apparent. 

A notable configuration with 4db/a 2 \/16.75 is shown in Fig. die), where the bonds that join 
the (lower corner of) honeycombs and the neighbors of the rebonded atoms are markedly stretched 
beyond their bulk value. The stretched bonds [indicated by arrows in Fig. Efc)] cause the surface 
energy to become very high, 100.39 meV/A 2 at the DFT level. Interestingly, the surface stress 
(and surface energy) is very efficiently lowered by adding dimers such that stretched bonds in 
Isjc) are broken, and the dimers "bridge" over as shown in Fig. E^a). The addition of one dimer 
per simulation cell results in the most stable Si(337) reconstruction that we found, with 7 ~ 89 
meV/A 2 . A structure similar to Eta) was previously proposed by Liu and coworkers ^3] to explain 
their HRTEM data. These authors also proposed a possible configuration for the other unit of (337) 
that is part of the large Si(5 5 12) unit cell, which has essentially the same bonding topology as|I](b). 
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Although no atomic scale calculations were performed in Ref. j!3j . the two (337) models of Liu et al. 
were in qualitative agreement with more recent DFT results by Jeong et al. |14| . as well as with STM 
and HRTEM images. When performing ab-initio relaxations, we found that the HOEP models that 
corresponded most closely to those in Ref. ^3] en d up relaxing into the configurations put forth by 
Jeong and coworkers [Tl]. The structures shown in Fig. 0a, b) allow for (337) nanofacets to intersect 
(113) facets without any bond breaking or rebonding. This absence of facet-edge rebonding, as well 
as the relatively low energy of Si(337), give rise to short-range attractive interactions between steps 
on miscut Si(113) 26 , which is consistent with the experimental observation of (337) step-bunched 
phases [TT] . 

5 Conclusions 

In summary, we have presented an investigation of the Si(337)-(2 x 1) reconstructions based on 
a genetic algorithm search for the most favorable atomic configuration |17l I21j . We have coupled 
the algorithm with the HOEP model [20] for atomic interactions in order to efficiently create a 
database of Si(337) models. Since no empirical interatomic potential can be expected to be fully 
transferable to arbitrary surface orientations, we have relaxed the database reconstructions using 
ab initio DFT calculations. The DFT relaxations give a more accurate idea about the relative 
stability of the reconstructions found, and help identify the likely lowest-energy candidates that 
are shown in Fig. 0(a) and (b). As it turns out, the two models coincide with those proposed 
by Jeong and coworkers (14] from STM experiments and density functional calculations. This 
finding lends strong support to the minimization approach taken here and shows that sampling 
the configuration space using the HOEP potential |20] is a versatile way to find good candidate 
reconstructions for high-index Si surfaces. We hope that the development of algorithms for surface 
structure determination ^JEl would also stimulate further research towards interatomic potentials 
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with improved transferability, since the use of such potentials will improve the quality of the model 
structures found during the global optimization: in turn, this would increase the range of materials 
problems where the search algorithms can be used to gain knowledge of optimal structure. 

The DFT surface energies of the (337) models shown in Fig. Efa,b) are relatively small, being 
about 1 meV/A 2 higher than the surface energy of Si(114) computed using the same computational 
parameters Section 2. Upon analyzing the results of the genetic algorithm optimization, we have 
seen that the (337) structures with low surface energies correspond to certain spatial arrangements 
of atomic-scale motifs (Fig. |2J), rather than to a minimum number of dangling bonds per unit 
cell. In experiments, particular conditions can readily appear such as to favor the stabilization of 
Si(337) over Si(5 5 12), as shown for instance in Refs. [31 H]. Recent works by Baski and coworkers 
|27| reveal a controlled way in which Si(337) facets can be formed -the deposition of gold. At Au 
coverages of 0.15 ML, the (5 5 12) orientation facets into Si(l 13) and Si(337), with the (337) surface 
unit cells being twice as large as the ones in Fig. ^ along the [776] direction. Since the Au coverages 
at which Si(337) emerges correspond to only a couple of Au atoms per each (337) unit, the models 
presented here can be used as building blocks for the structure of Si(337)-Au nanofacets observed 
in [27], as well as in other experimental situations. 
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Table I: Surface energies of selected Si(337)-(2 x 1) reconstructions, sorted by the number of atoms n 
in the periodic cell. The second column shows the number of dangling bonds per unit area, counted 
after relaxation with HOEP; the dangling bond density at the DFT level is shown in parentheses. 
Columns three and four list the surface energies given by the HOEP potential |2()j and by density 
functional calculations 22 . 



n 


Bond counting 


HOEP 


DFT 


Fig. 




(db/aVl6.75) 


(meV/A 2 ) 


(meV/A 2 ) 




266 


10 (10) 


87.37 


109.04 






10 (10) 


87.74 


97.22 


EDO 




6(4) 


87.79 


100.39 


He) 




10 (10) 


89.60 


102.35 






10 (10) 


89.79 


101.58 


He) 




(10) 




101.16 


BP) 




4(3) 


92.07 


94.47 


Ha) 




4(4) 


93.87 


101.43 


Btd) 


267 


8(8) 


92.35 


99.95 






8(8) 


92.37 


96.47 






10 (8) 


92.47 


107.72 






9 (8) 


99.11 


101.33 




268 


8(6) 


81.99 


93.13 


He) 




8(8) 


83.11 


95.17 






8(6) 


83.43 


89.61 


BM 




8(8) 


84.90 


94.52 


Hd) 




8(6) 


85.47 


90.19 






8(6) 


85.94 


89.12 


Ha) 


269 


4(4) 


89.18 


93.50 


Ha) 




6 (6) 


92.58 


99.21 


Hb) 
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Figure Captions. All figures COLOR ONLINE ONLY 

Fig. 1. Top and side view of the unreconstructed Si(337) surface. The rectangle represents the 
simulation cell that is subjected to shifted periodic boundary conditions and the solid arrows 
are the corresponding periodic vectors; four such periodic cells are shown in the figure in order to 
aid the eye with the shifted boundary conditions. The darker shade marks the undercoordinated 
surface atoms. 

Fig. 2. Structural features (top and side views) that can be present on low-energy Si(337) 
reconstructions: (a) dimers, (b) rebonded atoms, (c) tetramers, and (d) honeycombs. In each case, 
the atoms that make up the motif are shown in darker shade. 

Fig. 3. Models of Si(337) reconstructions with n = 266 atoms per (2 x 1) unit cell, after DFT 
relaxation (top and side views). The surface energy computed from first-principles is indicated for 
each structure, along with the corresponding value (in parentheses) determined with the HOEP 
interaction model [20]. The dark shade marks the undercoordinated atoms. The four-coordinated 
atoms that are exposed at the surface are shown in white. Apart from the relaxations, dimer tilting, 
and perhaps the relative phase of dimerization, the structure shown in panel (b) is the same as that 
proposed in Ref. |5|, and reconstructions in panels (e) and (f) are similar to the model in Ref. [S]. 

Fig. 4. Si(337) reconstructions with n = 268 atoms per (2 x 1) unit cell, after DFT relaxation 
(top and side views). The surface energy computed from first-principles is indicated for each struc- 
ture, along with the corresponding value (in parentheses) determined with the HOEP interaction 
model |2Uj . The dark shade marks the undercoordinated atoms, while the four-coordinated atoms 
that are exposed at the surface are shown in white. The lowest-energy Si(337)-(2 x 1) reconstruc- 
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tions are shown in panels (a) and (b). These structures differ in the position of the dimers in 
the unit cell: the dimers D can be part of a 7-member ring [side view, panel (a)] , or be part of a 
tetramer T [panel (b)]. The models shown in panels (a) and (b) are the same as those reported in 
Ref. ^1] for the two (337) units that are part of a (5 5 12) cell. 

Fig. 5. Models of Si(337) reconstructions with n = 269 atoms per (2 x 1) unit cell, after DFT 
relaxation (top and side views). The surface energy computed from first-principles is indicated for 
each structure, along with the corresponding value (in parentheses) determined with the HOEP in- 
teraction model [20]. The dark blue marks the undercoordinated atoms, while the four-coordinated 
atoms that are exposed at the surface are shown in white. The best n = 269 model [panel(a)], is 
stabilized by flat pentamers and subsurface interstitials. 



18 
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Figure 2: COLOR ONLINE ONLY. Chuang, Ciobanu, Wang, Ho 
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Figure 3: COLOR ONLINE ONLY. Chuang, Ciobanu, Wang, Ho 
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(a) y=89.12meV/A 2 (85.94) (b) y=89.61meV/A :: (83.43) 




3.13meV/A 2 (G1.99) (d) j =94.52 meV/A= (84.90) (c) y=9 



Figure 4: COLOR ONLINE ONLY. Chuang, Ciobanu, Wang, Ho 
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(a) y =93. 50 me V/ A 2 (89.18) (b) y =99.21 meV/A 2 {92.58} 
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